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Abstract 

We use the qualitative insight of a planar neuronal phase portrait to detect an excitability switch in arbitrary 
conductance-based models from a simple mathematical condition. The condition expresses a balance between 
ion channels that provide a negative feedback at resting potential (competitive channels) and those that provide 
a positive feedback at resting potential (cooperative channels). Geometrically, the condition imposes a trans- 
critical bifurcation that rules the switch of excitability through the variation of a single physiological parameter. 
The application of this algorithm to six different published conductance based models always finds the trans- 
critical bifurcation in a physiological parameter range, which suggests that the mathematical predictions have 
a physiological relevance and that a same regulatory mechanism is potentially involved in the excitability and 
signaling of many neurons. An experiment on midbrain dopaminergic neurons illustrates the predictive power 
of the balance condition by highlighting a role of cooperative L-type calcium channels in regulating the neuron 
firing pattern. 

Author summary 

Understanding the changing electrophysiological signatures of neurons in different physiological and pharmaco- 
logical conditions is a central focus of experimental electrophysiology because a key component of cell signaling 
in the nervous system. Computational modeling may assist experimentalists in this quest by identifying core 
mechanisms and suggesting pharmacological targets from a mathematical analysis of the model. But a suc- 
cessful interplay between experiments and mathematical predictions requires new analysis tools adapted to 
the complexity of high-dimensional computational models nowadays available. We use bifurcation theory to 
propose a mathematical condition that can detect an important switch of neuronal excitability in arbitrary 
conductance-based neuronal models and we illustrate its physiological relevance in six published state-of-the art 
models of different neurons. 

Introduction 

Detailed computational conductance-based models have long demonstrated their ability to faithfully reproduce 
the variety of electrophysiological signatures that can be recorded from a same neuron in varying physiological 
or pharmacological conditions. But the predictive value of a computational model is limited unless its analysis 
sheds light on the core mechanisms at play behind a computer simulation. Because conductance-based models 
arc nonlinear dynamical models, their analysis often requires a drastic reduction of dimension. The reduced 
model is amenable to the geometric methods of dynamical systems theory, but the mathematical insight is often 
gained at the expense of physiological interpretability; hence the need for methodological tools that can relate 
mathematical predictions of low-dimensional models to physiological predictions in detailed conductance based 
models. 
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In recent work [T] , we used phase plane analysis and dynamical bifurcation theory to characterize in reduced- 
order neurodynamics models a switch of excitability that is consistent with many physiological observations. 
More precisely, a transcritical bifurcation governed by a single parameter was shown to organize a switch 
from competitive excitability, extensively studied in most models inspired from the Hodgkin-Huxley model, to 
cooperative excitability, whose distinct electrophysiological signature is shared by most bursting neurons. 

The main contribution of the present paper is to show that this transcritical bifurcation, and the associated 
excitability switch, exist in a number of high- dimensional conductance-based models and that the resulting 
mathematical predictions have a strong physiological relevance. More precisely, we provide an algorithm to 
trace the transcritical bifurcation in arbitrary conductance-based models. Although purely mathematical in 
nature, the detection of the transcritical bifurcation is based on an ansatz that leads to a remarkably simple 
physiological interpretation: the switch of excitability is determined by a balance between competitive (those 
providing a negative feedback) and cooperative (those providing a positive feedback) ion channels at the resting 
potential. Because this simple tuning rule can take many different physiological forms, it is potentially shared 
by very different neurons. We apply the algorithm to detailed conductance-based models of six neurons known 
to exhibit drastic changes in their electrophysiological signatures depending on environmental conditions: the 
squid giant axon [2J, the dopaminergic neuron [3], the thalamic relay neuron 0], the thalamic reticular neuron 
[5], the aplysia R15 model [6], and the cerebellar granular cell [7]. In each case, the algorithm identifies a the 
transcritical bifurcation that occurs in a physiological range of parameters and its prediction correlates with 
experimental observations. 

In the case of dopaminergic neurons, our mathematical analysis suggests that the heterogeneity of L-type 
calcium channels across populations has little effect on the pacemaking but a strong effect on bursting. We 
validate this prediction via experimental results that are in strong agreement with the model prediction. 

After defining a novel classification of ion channels based on their competitive or cooperative nature, we 
briefly review the planar model presented in pQ and how its transcritical bifurcation qualitatively captures 
the switch between competitive and cooperative excitability. As a generalization of this low-dimensional case, 
we mathematically construct the same bifurcation in generic conductance based models and derive the balance 
condition determining the cooperative or competitive nature of the model. This construction and its electrophys- 
iological predictions are firstly illustrated on the squid giant axon. An algorithm for generic conductance-based 
models is subsequently derived and different models analysed. In particular, we revisit, using a TC neuron 
model as an insightful example, how neuronal signalling can rely on the dynamic regulation of the balancing be- 
tween competitive and cooperative ion channels. Finally, we illustrate our results in an experiment on midbrain 
dopaminergic neurons. 

Results 

Competitive and cooperative ion channels 

Conductance-based models of neurons describe the dynamic interaction between the membrane potential V and 
- possibly many - gating variables that control the ionic flow through the membrane. The gating of ion channels 
occurs on many different timescales. However, gating timescales can be grouped in three families, according to 
their influence on neuronal excitablity [8]: 

(i) Fast gating variables: These variables have a time-constant in the millisecond range. They generate 
the rapid regenerative upstroke of an action potential. Prominent representatives of this family are the 
activation gating variables of fast voltage-gated sodium channels (Nayl.l to Nay 1.9). 

(ii) (Slow) Recovery gating variables: These variables have a time constant 5 to 10 times larger than 
fast gating variables. They influence the spike initiation, downstroke, and the afterspike period. They are 
the key players of neuronal excitability. Prominent representatives of this family are the activation gating 
variables of delayed rectifier potassium channels (Kyl.l to Kyi. 3, Kyi. 5 to Kyi. 8, Ky2.1, Ky2.2, Ky3.1, 
Ky3.2, Ky7.1 to Ky7.5, KylO.l) and the activation gating variables of all calcium channels (Cay 1.x, 
Cay2.x, Cay3.x). 

(iii) (Ultra-Slow) Adaptation variables: These variables gate too slowly to be strongly activated by sin- 
gle action potentials. They modulate neuronal excitability only over periods of many action potentials. 
Prominent representative of this family are the inactivation gating variables of transient calcium channels 
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(Cay2.x, Cay3.x). Adaptation variables might also include non gating variables. For instance, the intra- 
cellular calcium concentration [Ca 2+ ]i„, which modulates the conductance of calcium-regulated channels. 

In view of their importance for neuronal excitability, we focus only on recovery gating variables to classify 
ion channels: when the recovery variable opposes a variation of the membrane potential by negative feedback, 
we term the associated channel a competitive ion channel. When the recovery variable instead enhances a 
voltage variation by positive feedback, the associated ion channel is termed cooperative. Ion channels that do 
not possess a recovery gating variable arc neither competitive nor cooperative and are called neutral . Neutral 
ion channels solely regulate the "quantity" of excitability without affecting its "quality". 

Table [I] shows a classification of many known ion channels according to this criterion. Not surprisingly, 
potassium channels are the main representatives of competitive ion channels. By increasing the total outward 
current, their activation induces a negative feedback on membrane voltage variations that is responsible for 
neuron repolarization. On the other hand, physiologically described calcium channels are all cooperative. 
Their activation induces an increase of the total post-spike inward current, in contrast to potassium channels. 
This is the source, for instance, of afterdepolarization potentials (ADP). Interestingly, sodium channels can be 
either competitive, cooperative, or neutral according to their fast transient, resurgent, or persistent behavior, 
respectively. 
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FA: fast activation, FI: fast inactivation, SA: slow activation, SI: slow inactivation, USA: utraslow activation, USI: 

ultraslow inactivation 



Table 1. Classification of ion channels according to their gating kinetics. Activation and 
inactivation variables are distributed in three groups: fast, (slow) recovery, and (ultra-slow) adaptation. 
Recovery variables are defined as competitive (resp. cooperative) if they induce a negative (resp. positive) 
feedback on membrane potential variations. An ion channel that posses a competitive variable is called 
"competitive channel", and similarly for cooperative channels. Channels that do not posses a recovery variable 
are called "neutral channels". This classification might change for a given channel for some channel subtypes. 

It is important to observe that the competitive (resp. cooperative) nature of channels is not solely linked 
to the outward (resp. inward) nature of the current. For instance, transient sodium channels are competi- 
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tive, because their recovery variable inactivates an inward current, inducing a negative-feedback on membrane 
potential variations. Similarly, potassium channels can be cooperative when their slow inactivation massively 
decreases the outward current, like in the case of A-type potassium channels. 

Although elementary, the classification above seems novel. It is motivated by the central message of this 
paper, that the balance between cooperative and competitive ion channels in a given neuron determines it 
neuronal excitability type. 



Competitive and cooperative excitability in planar models 



Planar models - that only consist of two state variables - have been instrumental to study excitability since 
the early days of neurodynamics 110] . Empirical planar reductions of conductance based models only retain 
the (fast) voltage variable V and one (slow) recovery gating variable n. Fast gating variables are set to steady- 
state (i.e. their fast variation is approximated as instantaneous), adaptation variables are treated as slowly 
varying parameters, and the sole recovery variable aggregates all recovery variables, expressing each of them as 
a (curve-fitted) state function of n. 

Based on such empirical reduction of the Hodgkin-Huxley model augmented with a calcium channel 
our recent study pQ explores the neuronal excitability of the planar model 

V = V - ^ - n 2 + I app (la) 

n = e(n o(V r - V ) + n Q - n) (lb) 

whose phase portraits are reproduced in Fig. [I] for two distinct values of the parameter uq (an indirect image 
of the calcium conductance in the high-dimensional model). The parameter e > is an image of the time- 
scale separation between V and n. The function rioo(-) has the standard sigmoid shape of conductance-based 
models and Vo is the half-activation potential. The phase portrait in Fig. [I] left is the qualitative phase 
portrait of FitzHugh-Nagumo model. It is associated to a reversible and sudden switch from rest to firing 
and has extensively been studied, with finer distinctions between excitability types I (SNIC), II (Hopf), III (no 
bifurcation), depending on the mathematical nature of the underlying bifurcation. See [T^l Section 3.4.4] .[T5].|141 
Section 7.1.3] and reference therein. 

The phase portrait in Fig. [T] right is in sharp contrast in that the electrophysiological response to a cur- 
rent input exhibits spike latency, plateau oscillations, and after depolarization potential (ADP). This specific 
signature, experimentally observed in many families of neurons, is fundamentally associated to the bistability 
illustrated in the phase portrait: namely, the robust coexistence of two stable attractors (a hyperpolarized 
resting potential and a limit cycle of periodic action potentials) and a saddle-separatrix that sharply separates 
their basins of attraction. The time evolution shown in the top figure is a consequence of this phase portrait 
and cannot be observed in FitzHugh-Nagumo like phase portraits. In that sense, it is an electrophysiological 
signature of competitive excitability. A finer distinction between excitability types IV (saddle-homoclinic) and 
V (saddle-saddle) is discussed in [1] depending on the nature of the bifurcations that surround the bistable 
range. 

A simple mathematical distinction between the two phase portraits shown in Fig. [T] is drawn from the 
Jacobian linearization of the model at the stable resting point (V , n): 
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The product of the partial derivatives = — 2n£^^-{V — Vo) is negative on the left phase portrait (small 

variations of the recovery variable provide negative feedback on the voltage) whereas it is positive on the 
right (small variations of the recovery variable provide positive feedback on the voltage). This difference is 
schematized in the block diagrams of Fig. [2] The reader familiar with mathematical biology will recognize the 
positive feedback interconnection as the basic circuit motif for mult i- stationary regimes |15j . 

The planar model can smoothly switch from competitive to cooperative excitability, with a transition 
occurring for n — 0, or, in algebraic terms, 

dV drioo 

d^mr- (2) 
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A convenient way to algorithmically track this excitability switch is to use bifurcation analysis and to impose 
that the critical condition ([2| coincides with a bifurcation of the model, which imposes the additional algebraic 
condition 



Simultaneously imposing ([2| and Q implies 



det J = 0. (3) 



which, in geometrical terms, corresponds to the transcritical bifurcation obtained for I app = | and illustrated 
in Fig. [3] 

The theory of bifurcation unfolding is further exploited in pQ in order to classify all excitability types 
associated to the planar model in Fig. [T] This analysis results in five different types of excitability obtained by 
varying the two parameters (Vo, no) around the singular phase portrait of Fig. |3j center. The parameter no acts 
in particular as the sole regulator of the balance between cooperative and competitive excitability by shifting 
the n-nullcline up and down: a positive n corresponds to a phase portrait as in Fig. [3] left, whereas the phase 
portrait of Fig. [3] right is obtained for sufficiently negative tiq. The associated parameter chart is reproduced 
in Figure [4j It contains the five types of excitability discussed above. 

Our main contribution in the present paper is to show that the diagram in Fig. [4] is not an artifact of planar 
reduction but captures excitability transitions that can algorithmically be tracked in conductance-based models 
of arbitrary dimension. 

Competitive and cooperative excitability in conductance based models 

We start by grouping gating variables of a given conductance-based model according to their time scales. The 
family Qp — {rriNa.f, niNa,p, tuk,a, ■ ■ ■} collects fast gating variables. The gating variable € [0, 1] denotes a 
generic member of this family. Similarly, the family Qr — {^./Va,/? tfik.dr? wiCa,L: • • 

.} collects recovery gating 

variables x r , whereas Qa = {hca.T, ^JVa,_R, "iif,M, • ■ ■} collects adaptation variables x a . For a given ion channel 
type i, the standard notation (resp. hi) is adopted for the activation (resp. inactivation) gating variable of 
the associated ionic current With these notations, a general neuron conductance-based model reads 

Cjr = -Y,9i< i h\ i {V-E i )+I avp , (5a) 

i 

T xf (V)xf= (xUV)-xf), (5b) 

r x r(V)x r = (x^V) - x r ), (5c) 

r x »(V)x a = (x a x (V)-x a ), (5d) 

where the sum in ([5p,) is over all ion channels in the model, and (|5p) , ^5p) , (|5p) hold for all the associated 
fast, recovery, and adaptation variables, respectively. The activation (resp. inactivation) functions x^, x 7 ^ are 
strictly monotone increasing (resp. decreasing) sigmoids. In the forthcoming analysis, all adaptation variables 
are treated as constant parameter, that is their slow evolution is neglected. 

We will detect a switch from competitive to cooperative excitability by mimicking the two-dimensional 
algorithm of the previous section. We first impose the bifurcation condition det J — 0, where J denotes the 
Jacobian of the subsystem (|^),([5b),([5p). The algebraic condition writes 

detJ / + V^=0 (6) 
1 Z-^dx r dV v ' 

where Jf denotes the Jacobian of the fast subsystem (|5£i) , ([JjJd) and the sum is over all recovery variables. The 
particular form of equation ([6]) is a direct consequence of the specific structure of conductance-based models, 
that is, parallel interconnection of two-dimensional feedback loops involving the voltage dynamics (|5^i) and one 
of the gating variable dynamics ((Hb) , (|5j;) . 

As for the planar model ([T]), we track the switch between competitive and cooperative excitability by imposing 
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the high-dimensional equivalent of the balance condition We therefore look for solutions of Q satisfying 
the ansatz 



E 



dV 3a£, = 
dx r dV 



(7) 



The two conditions now imply 



det Jf = 0. 



(8) 



We show in Supplementary Material that the corresponding bifurcation is necessarily transcritical. 

The singularity ([8| is the high-dimensional counterpart of the y-nullcline self-intersection in the planar 
model. It reflects the geometric nature of the transcritical bifurcation, that is, a robust geometrical object 
that exists independently of the timescale separation and persists in the singular limit of an infinite timescale 
separation, regardless of the system dimension. Our ansatz makes the proposed analysis completely robust 
against the model time constants. The time constants are only used to classify the gating variables in the three 
physiological groups. 

We split Gr in the two subfamilies Gr+ , which contains cooperative recovery gating variables x r+ , and Gr- , 
which contains competitive recovery gating variables x r ~ . The balance condition Q is then rewritten as 
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Cooperative gates 



Competitive gates 



to express a balance between competitive and cooperative ion channels. It is the high-dimensional counterpart 
of ([2]) and it provides a rigorous high-dimensional generalization of competitive and cooperative excitability: 



Competitive excitability 
The model is said to be competitive at steady state if: 



dv dx r + 

Z-^dx r + dV jL^i 
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Cooperative excitability 
The model is said to be cooperative at steady state if: 
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The insight provided by the planar model of the previous section predicts that the switch of excitability 
detected by the balance equation Q will lead to the accompanying distinct electrophysiological signatures 
of Fig. ' 
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Tracking excitability switches in the squid giant axon 

The Hodgkin-Huxley (HH) model [5] provides a non-physiological, but historical and experimentally verified 
tutorial for tracking a switch of excitability in conductance based models. The model reads 



CV = -g K n 4 (V-V K ) 
T m (V)m = (m oc (V r ) - m), 
r h (V)h={h oc {V)-h), 
T n {V)n = (noo(F) - n), 



g Na m 3 h(V -V Na )-gi(V-V l ) + I a 



(10a) 
(10b) 
(10c) 
(lOd) 



where m is the fast sodium channel activation while the sodium channel inactivation h and the potassium 
channel activation n are the (slow) recovery gating variables. We set all time constants to one, because this 
simplification has no effects on the algebraic conditions ^ and (JsJ) . The Jacobian of (10 1 reads 
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(11) 



The upper-left block is the Jacobian of the fast (V, m) subsystem. Imposing the singularity condition ^ yields 



dV 
dV 
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TC 



dm dV 



= 0, 



TC 



while the balance equation ^ reads 



dV drioo 



TC 



dV dh^ 
dh dV 



(12a) 



(13a) 



TC 



Note that (12) and (13) imply the bifurcation condition det J = in (11) 



At first sight, the balance condition ( 13 ) cannot be satisfied because both sodium and potassium channels are 
competitive channels according to their corresponding kinetics in the model, and in agreement with our proposed 
classification. This is consistent with the fact that the excitability of the HH model is always competitive in 
physiological conditions. 

However, it was long recognized [TB] that potassium channels can generate an inward current at steady- 
state if the extracellular K + concentration is sufficiently large. Indeed, any change in extracellular potassium 
concentration induces a change in the potassium reversal potential, as expressed by the Nernst equation. This 
suggests to use the potassium reversal potential Vk as a bifurcation parameter in HH model in order to satisfy 
the balance equation 



! >0 
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dV dn c 



dn dV 



dV dh. 



TC 



dh dV 



(14) 



TC 



n: Cooperative gating 



h: Competitive gating 



where potassium now acts as a cooperative gating variable provided that Vk > V$s- Physiologically, condi- 



tion ( 14 ) imposes that the potassium Nernst potential is large enough for the cooperativity of the potassium 



activation to balance the competitive effects of the sodium current inactivation. 



The two conditions ( 12 ) and ( 14 ) can numerically be solved to determine the critical values V K and V c 

which gives 



The 



value of the applied current at the transcritical bifurcation is then determined from ( 10 1 

A AV c )(V c -V K )+g Na m 



3 {V^h^V^iV 



V Na )+ 9l {V c -Vi). 
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The numerical bifurcation diagram in Fig. [j^V confirms the transcritical bifurcation at V K . That bifurcation 
diagram is drawn by varying Vk together with applied current I app , following the affine reparametrization 
described in Supplementary Material. More precisely, 

I apP (V K ) = I c app - g K nUV c )(VK - V K ). 

Mathematically, this reparametrization imposes one of the defining conditions of the transcritical bifurcation. 
Physiologically, its effect is to keep the net current constant at steady-state V c (i.e. Iion(Vic) + Ih PP (Vk) = 0): 
as Vk is varied, the observed switch in the excitability type does not rely on changes in the net current across 
the membrane, but solely on changes in its dynamical properties. 

The bifurcation diagram in Fig. [5]A_ provides informations on the model excitability also far from the 
transcritical values. For highly hyperpolarized Vk, the model is purely competitive and exhibits the typical 
excitable behavior of the original Hodgkin-Huxley model, which corresponds to Type II excitability (see Figure 
[4] and pfl Figure 8]). As Vk is increased, a stable cooperative steady-state is born in a saddle-node bifurcation. 
At this transition, the system switches to Type V excitability, which is cooperative. Short current pulses let 
the system switch between the depolarized competitive steady state and the hyperpolarized cooperative stable 
steady state (Fig. [5j3, middle). The bifurcation diagram and the phase-portrait carry the signature of Type V 
excitable models (Fig. |5p,D, middle) [H Section 5]. Finally, further increase of Vk let the competitive steady 
state exchanges its stability with a (cooperative) saddle at the transcritical bifurcation and the system switches 
to Type IV excitability. The cooperative steady state coexists in this case with the spiking limit cycle attractor. 
Current pulses switch the systems between the two attractors (Fig. |5p, right), and the bifurcation diagram and 
phase portrait carry the signature of Type IV excitable models (Fig. [5p,D, right) [TJ Section 5]. 

The same qualitative excitability switch was described by Rinzel in [10] . who linked the appearance of a 
bistable behavior to the inward nature of potassium current at steady-state for sufficiently depolarized Vk- In 
vitro recordings of the squid giant axon with isotonic extracellular K + concentration show the same transition 
|16| . Our interpretation of this phenomenon in terms of competitive and cooperative excitability is as follows: 
as the potassium Nernst potential is increased, the potassium current activation gating variable becomes less 
and less competitive and eventually becomes cooperative. Further increases of Vk causes its cooperativity to 
dominate the competitivity of the sodium current inactivation, providing an overall positive feedback on the 
fast subsystem composed by the membrane potential and the sodium current activation, leading to cooperative 
excitability, as sketched in Fig. [6] 

Tracking excitability switches in conductance-based models 

The mathematical analysis of the previous sections follows an algorithm that allows to detect a transcritical 
bifurcation in generic conductance based models of arbitrary dimension and to track associated excitability 
switches. The steps of the algorithm are summarized in Table [2] 

We now apply this algorithm to a number of published conductance-based models and show that all these 
models can switch between competitive and cooperative excitability through a transcritical bifurcation, as 
sketched in Fig. [7] Figure [7] indicates two qualitatively distinct routes from competitive to cooperative ex- 
citability: the II-V-IV route just described with Hodgkin-Huxley model (Fig. [5j\.) and the II-I-IV route that 
will be illustrated on the dopaminergic neuron model. 

Dopaminergic (DA) neuron model 

Model equations and parameters are taken fom This is the model that originally motivated [IT] . 

(i) Classification of gating variables as fast (Gf), recovery (Gr), and adaptation (Ga) variables 
The model includes fast sodium channels (Ino,j), delayed-rectifier potassium channels (Ik, dr.), L-typc 
calcium channels (Iccl,l), small conductance calcium-activated potassium (SK) channels (Ik.Co), and 
calcium pumps (Ica,jmm P )- We classify model variables as follows 

• Gf = {m Na j} 

• Or - = {hffa.f, m K,DR\ and Gr,+ = {m Ca x} 

• Ga = {[Ca 2 +] m } 
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Algorithm for the detection of a transcritical bifurcation in generic conductance-based models 
and computation of the excitability switch bifurcation diagram. 

(i) Classification of gating variables as fast (Gf), recovery (Gr), and adaptation (Ga) variables 
(i-a) Following Tab. [TJ group gating variables in the three groups Gf, Gr, and Ga- 
(i-b) Split Gr, in cooperative Gr,+ and competitive Gr,- recovery gating variables, 
(i-c) If adaptation variables are present, set them to constant physiologically relevant values. 



(ii) Balance equation and choice of the bifurcation parameter 

(ii-a) Select a cooperative ionic current I coop and the associated cooperative recovery gating 
variable x coop . 

(ii-b) Isolate the contribution of I coop to excitability by setting to zero the conductances of all 

other cooperative currents, if any. 
(ii-c) Write the balance equation 



dV dx c °°P 



Q x coop gy 



dV dxlr 



TC 

x 



dx r dV 



(b.eq.) 

TC 



(ii-d) If I coop has an adaptation variable x a , pick it as the bifurcation parameter A regulating 
the left hand side of (b.eq.), that is A = x a . 
If I C oop has no adaptation variable, pick A = g COO p- 
(iii) Singularity condition and fixed point equation 

(iii-a) Solve (b.eq.) together with the singularity condition dsj), in V and A. 

For numerical implementation, note that the left hand side of QSp can be written as 

d 

I l 

cf=xl (V) / 

(iii-b) Plug the computed values V c and A c into the fixed point equation V\tc = to 
compute the value of the applied current at the transcritical bifurcation (/ c ). 



det Jf = — I V 



(iv) Tracking of excitability switches 

Change the applied current according to the equation 

dV 

(A-A c ), 



!app - I q x 



TC 

and compute the model bifurcation diagram with V as the variable and A as the bifurcation parameter. 



Table 2 
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In order to unfold excitability switches, SK channel density is set to zero, since SK channels drastically 
attenuate DA neuron excitability by activating a strong calcium regulated potassium current jTTl [TBI E] ■ 
The intracellular calcium concentration is fixed at [Ca 2+ ]i n — 300 nM. 

(ii) Balance equation and choice of the bifurcation parameter 

The only source of cooperative excitability is provided by L-type calcium channels. The balance equation 
reads 



dV 9TO C a,L,oo 



dm C a,_ 



dV 



dV dh 



Na,f,<x 



dV dm Ki DR,oo 



TC 



TC 



\dh Na j dV dm K ^ DR dV J 

We use the L-type calcium channel density gca.L as the bifurcation parameter (i.e. A = gca,h) 
Solving Steps (iii) and (iv) of the algorithm above gives the following resulttj^J 



(15) 



v c 


9Ca,L 


r 


-64.9167 mV 


1.9473 mS/cm 2 


9.5436 ^A/cm 2 



Note that all values are within the physiological range of the original model. In particular, the critical value 
9ca l i s roughly 1.5 times smaller than the original value one, which is consistent with the observation that the 
original model exhibits cooperative excitability during SK channel blockade. 

The resulting bifurcation diagram is drawn in Fig. [HJ\. In addition to confirming the existence of a transcritical 
bifurcation for the computed values, it reveals the excitability switches induced by changes in L-type calcium 
channel density in this model: in the absence of L-type calcium channels, the model exhibits Type II excitability. 
As gca.L increases, a saddle point and an unstable node emerge at a saddle-node bifurcation: the model switches 
to type I excitability. Further increase of gca,L causes a transcritical bifurcation, where the stable point and 
the saddle exchange their stability. At this point, the stable steady-state becomes cooperative, and the model 
switches to cooperative Type IV excitability. 

These excitability switches induce the predicted changes in the electrophysiological signatures, as illustrated 
in Fig. [8j3. Whereas the DA neuron model instantaneously reacts to a step input of depolarizing current for 
gCa,L < (ica l ' ^ exhibits electrophysiological signature of cooperative excitability such as spike latency, plateau 
oscillations and ADP as soon as gca.L becomes higher than g c Ca L . In addition, the model becomes strongly 
bistable. 

Thalamic relay (RE) neuron model 

Model equations and parameters are taken from 0]. 

(i) Classification of gating variables as fast (Gf), recovery (Gr), and adaptation (Ga) variables 
The model includes fast sodium channels -T/Va,/, delayed- rectifier potassium channels Ik,dr and T-type 
calcium channels Ica,T- We classify model variables as follows 

• Gf = {niNaj} 

• Gr - = {h Na j, m K , DR } and Gr,+ = {m Ca ^ T } 

• Ga = {hca,T} 

(ii) Balance equation and choice of the bifurcation parameter 

The only source of cooperative excitability is provided by T-type calcium channels. The associated balance 
equation reads 

dV dm C a,T, 
dm Ca ,T dV 

T-type calcium channels are dynamically regulated by a (slow) voltage-gated inactivation hc a ,T- We use 
this variable as the bifurcation parameter (i.e. A = /ica.r)- 

1 In this and the subsequent tables, the reported critical values are not exact, but rounded to the last significant digit. 
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Solving Steps (iii) and (iv) of the algorithm above gives the following results: 









v c 


' L Ca,T 


r 


-61.5561 mV 


0.0041 


0.9074 nA/cm 2 



Again, these values are in a physiological range, and the electrophysiological signatures are consistent with the 
excitability switches. See Fig. [9] 

Thalamic reticular (RT) neuron model 

Model equations and parameters are taken from [5], maximal conductances are adapted as in |llj . 

(i) Classification of gating variables as fast (Gf)i recovery (Gr), and adaptation (Ga) variables 
The model includes fast sodium channels Ino,j, delayed- rectifier potassium channels Ik.dr and T-type 
calcium channels Ica,T- We classify model variables as follows 

• Gf = {m N a,f} 

• Gr - = {hffa.f, m K ,DR} and Gr,+ = {m Ca ,T} 

• GA= {hca,T} 

(ii) Balance equation and choice of the bifurcation parameter 

The only source of cooperative excitability is provided T-type calcium channels. The associated balance 
equation has the same structure as for the thalamic relay neuron model considered above. Along the same 
line, we choose the T-type calcium channel inactivation hc a ,T as the bifurcation parameter (i.e. A = hr). 

Solving Steps (iii) and (iv) of the algorithm above gives the following results: 









v c 


,l Ca.T 


r 


-48.8063 mV 


0.1780 


-0.82484 nA 



As in the case of the thalamic relay neuron model, the T-type calcium channel inactivation generates a dynamical 
switch between competitive and cooperative excitability, significantly affecting neuron response to external 
inputs (Fig. [9j). 

Aplysia R15 neuron model 

Model equations and parameters are taken from [BJ. 

(i) Classification of gating variables as fast (Gf)i recovery (Gr), and adaptation (Ga) variables 
The model includes fast sodium channels InhJ, delayed-rectifier potassium channels Ik,dRi slow L-typc 
calcium channels Ica,L and calcium-activated potassium channels Ik,Co,- We classify model variables as 
follows 

• Qf = {m N a,f} 

• Gr - = {hNaj, m K ,DR} and Gr,+ = {m Ca x} 

• Ga = {[Ca 2 +] m } 

The intracellular calcium conductance is fixed at [Ca 2+ ]i n = 0.09?iAf . 

(ii) Balance equation and choice of the bifurcation parameter 

As in the case of DA neurons, the source of cooperative excitability is provided by L-type calcium channels, 
and we take their maximal conductance gca.L as the bifurcation parameter. The associated balance 
equation has the same structure as for the DA neuron model considered above. 
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Solving Steps (iii) and (iv) of the algorithm above gives the following results: 



v c 


9(Ja,L 


r 


-48.0516 mV 


5.4054 10~ 5 mS/cm 2 


-0.0318 fiA/cm 2 



Comparing the critical value g c Ca L with the original value gca,L = 4 10~ 3 mS/cm 2 shows that the bursting 
model proposed in [5] exhibits strong cooperative excitability. Switches of electrophysiological signatures are 
illustrated in Fig. [9j 

Cerebellar granular cell ( GC) model 

Model equations and parameters are taken from [7]. 

(i) Classification of gating variables as fast (Gf), recovery (Gr), and adaptation (Ga) variables 
The model includes the following ion channels: fast (iNa.f), persistent (Inci,p) and resurgent sodium 
channels (Ino,,r) ', N-type calcium channels (Ico,,n) ; delayed rectifier (Ik.dr), A-type (Ik, a), inward 
rectifier (Ik,ir)i calcium activated (Ik, Co) an d slow potassium channels (Ik,sIow)- We classify model 
variables as follows 



• Gf = {m-Naj, rriNa,p} 

• Gr - = {hffaj, m K ,DR, m K ,A, m K ,iR\ and Gr.+ 

• Ga = {hNa,R, hca,N, h>K,A, mK,slow} 



{mNa,R, m C a,N} 



We set the persistent and calcium-activated currents to zero (these two channels do not impact excitability 
type as anticipated by our classification and shown by D'Angelo and colleagues [7]). The inactivation of 
the A-type potassium current is fixed at Kk,a = 0.02 and the activation of the slow potassium current is 
fixed at rriK,siow = 0.13. 

(ii) Balance equation and choice of the bifurcation parameter 

The neuron model possesses two sources of cooperative excitability: resurgent sodium channels and N-type 
calcium channels. We apply our algorithm considering each case at a time, setting the other cooperative 
source to zero. As in the case of T-type calcium channels mentioned above, these two channels possess 
an inactivation gate, which is used as the bifurcation parameter. The associated balance equations read, 
respectively: 



dV drriNa 



R,oo 



9m Na>R dV 

dV dm C a,N,oo 
dm Ca ,N 9V 



dV dh Na 



f,°° 



TC 



TC 



dh NaJ dV 

dV dh NaJ/ , 
dh NaJ dV 



+ 



E 

DR, A, 

E 



dV drriK.y, 
dm,K. v dV 



y=DR, A, IR 



y=DR, A, IR 



dV dm K , y ,oo 
dm,K, y dV 



TC 



TC 



Solving Steps (iii) and (iv) of the algorithm above gives the following results, respectively: 



v c 


a Ca,N 


r 


V c 


h c 

n Na,R 


r 


-41.5447 mV 


0.0318 


5.2508 fiA/cm 2 


-41.5252 mV 


0.5066 


5.1549 //A/cm 2 



These results show that both channel densities are sufficiently high to allow dynamical excitability switches, 
which are due to dynamical channel (de)inactivations. Note that the N-type calcium channels contribute much 
more to cooperative excitability than resurgent sodium channels in this model (h c Ca N « h c Na R ). The effect 
of N-type channel regulation on GC excitability are shown in Fig. [9j A similar switch is observed for resurgent 
sodium channels 
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Neuronal signaling through dynamic balancing of cooperative and 
competitive ion channels 

In several of the examples analyzed in the previous section, the bifurcation parameter is a voltage-dependent 
(ultra-slow) adaptation gating variable of a cooperative ion channel: the T-type calcium channel inactivation 
hca,T, the N-type calcium channel inactivation hc a ,N, or the resurgent sodium channel inactivation hjq a ,R- The 
physiological implication is that, in all these examples, any regulation of the resting potential acts as a regulator 
of the balancing between cooperative and competitive ion channels and therefore as a regulator of excitability. 

We illustrate this regulation mechanism on reticular and relay cells of the thalamus, whose firing patterns 
have been linked to awakeness and sleep [2U] . 

Fig. [10] A. illustrates that, as predicted by our bifurcation analysis, the excitability of RT cells is cooperative 
at hyperpolarized resting potential of —95 mV (hca.T ~ 1 > h,Q a T ) whereas it is competitive at a normal 
resting potential of —57 mV (hc a ,T <C h c Ca T ). 

The amount of cooperative channels that can be activated therefore depends on the resting potential. The 
latter can be regulated by GABAg metabotropic receptors, for instance. As a consequence, the model predicts 
two possible excitability patterns in response to external sensory stimuli: 

(i) "Linear integration" (competitive) - With a depolarized resting potential, the competitive (Type I) ex- 
citability of the neuron causes a linear integration of external stimuli, faithfully transmitted to the cortex. 
In the case of the lateral geniculate nucleus, this firing mode has been suggested to support a faithful 
reconstruction of the visual world |21j . 

(ii) "Wake-up call" (cooperative) - With a hyperpolarized resting potential, the cooperative (Type IV) ex- 
citability of the neuron is associated to a burst response to external stimuli. The signal-to-noise ratio is 
higher during bursting, which enhances its detectability. This maximization of stimulus detection has been 
suggested to detect sudden changes in the environment |21j . The bursting mode is transient because of the 



slow inactivation of calcium channels (Fig. 10 3), letting the model switch back to the linear integration 
competitive regime. 

This prediction is consistent with the hypotheses that the transient burst could serve as a "wake-up call", 
switching the neuron from a silent state to a responsive state where stimuli can be transmitted to the cortex. 
These two firing modes correspond to the "up-state" of the thalamus, where external stimulations are detected, 
integrated and transmitted to the cortex. Note that the "down-state" of the thalamus, observed during slow- 
wave sleep and absence seizures, is associated to spontaneous synchronized burst of TC neurons. Those cellular 
rhythms also exhibit cooperative excitability but the synchronization of bursts involves additional population 
mechanisms not discussed in the present paper. 

In conclusion, our analysis suggests that the complex firing patterns of the TC neuron owe much to a simple 
but regulated balance between competitive and cooperative channels. 



DA neuron signaling and cooperative excitability 

Our analysis in the previous sections predicts that the L-type calcium conductance gca.L is a tuner of cooperative 
excitability in the DA neuron model [3] . Specifically, our model predicts cooperative excitability in the absence 
of small conductance calcium-activated potassium (SK) channels, provided that gca.L is sufficiently large. 

Physiologically, the parameter g~ca,L is an image of the L-type calcium channel density. The density of 
these channels has been suggested to exhibit some heterogeneity among DA neuron population [231 1231 ISj ■ As 
a consequence, the type of excitability of DA neurons under the blockade of SK channels is expected to be 
heterogenous across a population of DA neurons. 

We tested this hypothesis on DA neuron extracellular recordings from the substantia nigra pars compacta 
(SNc) in the presence of synaptic blockers (10 mM CNQX, 1 mM MK801, 10 mM SR95531, 1 mM sulpiride 
and 1 mM CGP55845, which block AMPA, NMD A, GABA A , D2 and GABA S receptors, respectively) (Fig. 



111. Control experiments showed that application of the synaptic blockers alone induced a small increase in 
firing rate which was stable for at least one hour The application of Apamin 300 nM, a SK channel blocker, 
has mixed effects on DA neuron excitability: whereas it barely affect firing rate and pattern in some neurons 
(Fig. [TT]k, left), it strongly increases irregularities and/or induces spontaneous bursting in others (Fig. [TT]A, 
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right). The application of a L-type calcium channel blocker (Nifedipine 20 fiM) reversed the effect of SK channel 
blockade on DA neuron signaling, regularizing the firing pattern in all neurons (Fig. |11| 3). This suggests that 
L-type calcium channels play an important role in the SK channel blockade induced irregularities and bursting 
of DA cells in vitro. 

These results show that the L-type calcium channel density has a critical impact on the spiking pattern and 
therefore on the excitability type of DA neurons, although it barely affects the pacemaking activity. If the L- 
type calcium channel density is low, the total ionic current provides a negative feedback on membrane potential 
variations (Fig. [Tip, top), and the neuron fires in tonic single spike, a signature of competitive excitability. 
By contrast, if the L-type calcium channel density is high, the positive feedback provided by these channels 
is dominant (Fig. |llp , bottom), and the neuron exhibits burst firing, a signature of cooperative excitability. 
The experiment confirms that in this latter situation only the application of the L-type calcium channel blocker 
nifedipine can have a critical impact on the firing pattern by inhibiting the positive feedback provided by calcium 
channels, therefore generating a switch from cooperative to competitive excitability. 

The role of L-type calcium channels in DA firing patterns has been controversial to date [55J 123] , mainly 
because the heterogeneity of L-type channels does not seem to affect the pacemaking activity of DA neurons. 
We explained in previous work 3 that this robustness relies on the collaboration of L-type calcium channels 
and sodium channels in generating the pacemaking activity. In this sense, calcium and sodium channels are 
redundant. In contrast, our analysis suggests that calcium channels and sodium channels have a different, and 
in fact antagonistic, contributions to the quality of neuronal excitability, a difference that is unravelled only 
by the blockade of SK channels. SK channels activate indeed a strong calcium regulated after-spike potassium 
current that masks any excitability switch by regularising the spiking pattern to a robust pacemaking activity 
[T71 ITS1 ITU] . To sum up, the L-type calcium channels strongly affect excitability and bursting as cooperative 
channels, whereas they can be substituted by sodium channels for pacemaking. 

Blockade of L-type calcium channels has recently been suggested as a potential guard against neurodegen- 
eration [2H [22]: it would protect the neuron against a calcium inflow associated to pacemaking, a source of 
metabolic stress favoring neurodegeneration [551 HZ] ■ Results in this section suggest however that this blockade 
would at the same time significantly affects the neuron excitability even if it does not affect the pacemaking 
activity. 

Discussion 

A simple and robust balance equation identifies a transcritical bifurcation in arbi- 
trary conductance based models 

Motivated by a geometric analysis of a qualitative phase portrait, we have proposed an algorithm that can 
detect a transcritical bifurcation in arbitrary conductance based models. Owing to the special structure of 
such models, the algorithm leads to solving an algebraic equation of remarkable simplicity and physiological 
relevance: the type of excitability of a neuronal model is primarily determined by a balance between competitive 
and cooperative ion channels. The condition is also robust because the balance is independent of the detailed 
kinetics. It only relies on a classification of variables in three well separated time-scales, in full accordance with 
what is known on ion channels kinetics [8]. 

A same mathematical prediction applies to many distinct physiological observation 

Although purely mathematical in nature, the transcritical bifurcation has a remarkable predictive value in several 
published conductance based models. In each of the six analysed models, the proposed algorithm identifies a 
physiological parameter that acts as a tuner of neuronal excitability in a physiologically plausible range and 
in full agreement with existing experimental data. At the same time, the distinct nature of the regulating 
parameter, which can be either the maximal conductance or the inactivation gating variable of a cooperative 
ion channel depending on the neuron model, is associated to distinctly different regulation mechanisms. 

Reduced modelling should retain the balancing channel 

The classification of gating variables in three distinct time scale is an essential modelling step both for the 
proposed algorithm and for the reduction of full conductance-based models to low-dimensional models that can 
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be used in population studies [28]. When all recovery ion channels are properly identified, they can be aggregated 
in a single recovery variable to lead to a second order model of the type ([l]), where the single parameter hq 
captures the competitive or cooperative nature of the aggregated recovery variable. Further reduction to a 
one-dimensional hybrid model with reset is possible thanks to the time-scale sepration between the voltage V 
and the recovery variable n. This reduction is illustrated in [11] on the thalamic TC neuron and the reduced 
model remarkably retains the switch of excitability of its high-dimensional counterpart. In contrast, a reduced 
model will loose the switch of excitability of the full conductance-based model when a cooperative ion channel 
is treated as a fast gating variable. This is for instance the case in the STN neuron model of [35] . The switch of 
excitability is well identified in experimental data [30] but cannot be reproduced in the model because the only 
source of cooperativity, the activation of T-type calcium channels, is treated as a (instantaneous) fast variable. 

Neuronal excitability is regulated 

In each of the analysed conductance-based models, the balance equation responsible for the switch of ex- 
citability is satisfied in a physiological range of parameters. This observation supports the hypothesis that 
neuronal excitability is highly regulated by molecular mechanisms. Our experiment on the dopaminergic neu- 
ron illustrates how the density of L-type calcium channels is a key regulator of excitablity. This density could 
be adapted, for instance, by phosphorylation/dephosphorylation mechanisms. From a physiological point of 
view, a potential regulation mechanisms of DA neuron bursting through L-type calcium channel phosphoryla- 
tion/dephosphorylation might rely on metabotropic receptors. 

Materials and Methods 

Numerical analysis 

Numerical temporal traces of the different neuron models were run in MATLALrl The phase portraits in Figures 
were hand-drawn using the Open Source vector graphics editor Inkscaper] The phase portraits in Figure 
[5] were numerically drawn with MATLAB and modified with Inkscape. The bifurcation diagrams in Figures [5] 
and [8] were drawn by implementing the algorithm of Table [2] in MATLAB. 

In vitro extracellular recordings 
Housing 

Adult (200-250 g) male Wistar rats were housed in groups of three or four, supplied with food and water ad 
libitum, and maintained on a 12 hour light/dark cycle. 

Tissue preparation 

Rats were anesthetized with chloral hydrate (400 mg/kg, i.p) and decapitated. The brain was extracted out 
of the skull in less than 1 min and cooled in ice-cold (2 — 4°C) artificial cerebrospinal fluid (ACSF) of the 
following composition (in mM): NaCl, 130; KCl, 5; NaH 2 P0 4 , 1.25; MgS0 4 , 1.25; CaCl 2 ,2; D - glucose, 10; 
NaHC0 3 ,2A; saturated with 95%0 2 and 5%C0 2 (pH 7.4). 

Horizontal slices (400 /im) containing the midbrain were cut in cold ACSF in a vibratome and transferred 
to a small beaker containing ACSF at 32°C. After 10-30 min slices were transferred to the recording chamber 
(volume 0.5 ml) in which they were completely submerged and superfused by ACSF at a rate of 2.5 ± 0.5 
ml/min (T°: 34.5 ± 0.5°C). The tissue was held in position by small platinum weights. 

Action potential recordings and identification of neurons 

Action potentials (amplitude: 200 to 1000 fiV) were recorded extracellularly using borosilicate capillaries (model 
1403516, Hilgcnberg, Malsfeld, Germany) which had been pulled on a Narishige (PE-2) vertical puller. The tip 
resistance was 5-20 MCI. They were passed through an impedance adapter and amplified one thousand times 
with a home made amplifier. They were displayed on an oscilloscope and fed to an analog digital interface 

1 http : //www . mathworks . com 
3 http : //inkscape . org 
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(CED 1401) connected to a computer. Data were collected with the use of the ?Spike 2? software (version 4.23, 
Cambridge Electronic Design, Cambridge, UK). Use of an appropriate template ensured that the activity of 
only one cell was recorded. 

Electrophysiological criteria were used in order to identify DA neurons. In vitro, these neurons exhibit a 
very regular firing pattern and long (> 2.5 ms), triphasic spikes, often displaying a prominent notch in the initial 
positive rising phase. They have a slow firing rate comprised between 0.5 and 5 Hz. A large number of previous 
experiments in our laboratory have shown that the firing of these cells is inhibited by nanomolar concentrations 
of D2 agonists such as quinpirolc and BHT920. This was verified again in this study in the absence of synaptic 
blockers. Neurons meeting these criteria were selected for recording. 

Only one cell was studied per slice. All recordings were made in the substantia nigra pars compacta area. 
Drugs were applied by superfusion at known concentrations using three-way taps. Each concentration was 
applied for at least 10 min to ensure that the drug concentration reached equilibrium in the tissue. Drug 
effects were quantified by calculating the mean firing frequency over the 5 last minutes of the control period 
and over the last minute during which the drug was superfused. Synaptic blockers were superfused throughout 
all experiments unless stated otherwise. Neurons whose firing rate varied by more than 5% (average of 1 min) 
during the control period were discarded. 

Experimental protocol 

The activity of DA neurons was recorded in control conditions for 5 minutes. Synaptic blockers (10 \iM CNQX, 
1 fj,M MK801, 10 [iM SR95531, 1 \xM sulpiride and 1 fiM CGP55845 to block AMPA, NMDA, GABA A , D2 
and GABAs receptors, respectively) were then superfused for the remainder of the experiment. After at least 10 
minutes, SK channels were blocked by bath application of 300 nM apamin for the remaining of the experiment. 
At least 10 minutes after the beginning of SK channel blockade, L-type calcium channels were blocked by bath 
application of 20 \iM nifedipine for 10 minutes as well. 

Drugs 

The sources of the drugs used were as follows. CGP55845, CNQX, MK801, TTX and Nifedipine were obtained 
from Tocris Cookson (Bristol, UK). SR95531 and apamin was obtained from Sigma (St Louis, MO, USA). 
Sulpiride was a gift from Sanofi-Aventis. 
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Figure 1. The bottom figures illustrate the typical phase portrait of competitive (left) and cooperative (right) 
excitability. The dark (resp. light) blue circle denotes a stable competitive (resp. cooperative) steady state 
(V,n). The full (resp. dashed) line is the voltage (resp. recovery variable) nullcline. The saddle point in the 
right phase portrait is represented as the cross and its separatrix as the green oriented line. The stable limit 
cycle surrounding the unstable fixed point (represented as a circle) is represented by the blue oriented line. 
The top figures illustrate the typical accompanying electrophysiological responses to step variations of current. 
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Figure 2. Block diagrams representation of competitive and cooperative excitability in planar models. 
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Figure 3. A continuous deformation from the competitive phase portrait of Fig. [lj left to the cooperative 
phase portrait of Fig. [I] right involving a transcritical bifurcation determined by the algebraic conditions ([2| 
and Q. The dark blue circle represents a competitive stable steady-state, the light blue circle a cooperative 
stable steady-state, and the half-filled circle represents the transcritical bifurcation which separates the 
competitive and cooperative regimes. 




Figure 4. Excitability types in model SN denotes the saddle-node bifurcation, TC the transcritical 
bifurcation. ■: Pitchfork bifurcation organizing center. Varying n Q and Vq the model switches between 
excitability types. 
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Figure 5. Variations of the potassium reversal potential Vk induce excitability switches in the 
Hodgkin-Huxley model. A. Bifurcation diagram of the HH model with Vk as the bifurcation parameter. 
TC denotes a transcritical bifurcation, SN a saddle-node bifurcation, HB a Hopf bifurcation. Branches of 
stable fixed points are represented as solid line, whereas branches of saddle points and unstable points as 
dashed lines. B. Electrophysiological responses of the model for three different values of Vk, corresponding to 
three different excitability types (Type I, V and IV, from left to right). C. Bifurcation diagrams with the 
applied current as the bifurcation parameter for the same three values of Vk as in B. Black (resp. blue) full 
lines represent branches of stable steady-states (resp. limit cycles), black dashed lines branches of saddle and 
unstable steady-states. Branches of unstable limit cycle are drawn as dashed blue lines. HB denotes a Hopf 
bifurcation, SN a saddle- node bifurcation, and SH a saddle- homoclinic bifurcation. D. Phase portraits of 
reduced HH model proposed by Rinzel in [TU] for the same three values of Vk as in B,C. Blue full lines denote 
the y-nullclines and black full lines the w-nullclines, where w denotes the recovery variable of the reduced 
model. Filled circles denote stable steady-states, crosses saddle points, and circles unstable steady-states. 
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Figure 6. Representation of HH model excitability as feedback interconnections. The membrane 
potential V and the sodium channel activation m are the fast variables, the potassium channel activation n 
and the sodium channel inactivation h are the slow recovery variables. When Vk is below (resp. above) the 
critical value V£, recovery variables induce an overall negative (resp. positive) feedback on membrane 
potential variations, resulting in competitive (resp. cooperative) excitability. 




Figure 7. Modifications in the balance between competitive and cooperative channels induce 
excitability switches in conductance-based models. The figure sketches excitability switches of the 
Hodgkin-Huxley (HH) model [2J, Aplysia's R15 neuron (R15) model [5J, a dopaminergic (DA) neuron model 
[3], thalamic reticular (RT) and relay (RE) neuron models [HE], an d a cerebral granule cell (GC) model [7j on 
the excitability parameter map computed for the two-dimensional model of pQ. All these conductance-based 
models can switch between competitive and cooperative excitability through the physiologically relevant 
regulation of specific ion channels. 
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Figure 8. Variations of L-type calcium channel density gca,L induces excitability switches in a 
model of DA neurons |3j. A. Bifurcation diagram of the model with g~ca,L as the bifurcation parameter. 
TC denotes a transcritical bifurcation, SN a saddle-node bifurcation. Branches of stable fixed points are 
represented as solid line, branches of saddle points and unstable points as dashed lines. B. 
Electrophysiological responses of the model to step inputs of excitatory/inhibitory current (the intracellular 
calcium concentration is fixed at [Ca 2+ ]i n = 300 nM, which is within the physiological range). For gca,L lower 
(resp. higher) than the critical value g^ a L , the model exhibits typical electrophysiological signature of 
competitive (resp. cooperative) excitability. The low gca,L configuration corresponds to gca.L — lmS/cm 2 , 
whereas the high g~ca,L configuration corresponds to g~ca,L = 3mS / cm 2 . 
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Figure 9. The same mathematical bifurcation in different conductance-base models causes the 
same switch in electrophysiological signatures. The figure shows the electrophysiological responses of 
various conductance-based models to step inputs of excitatory/inhibitory current when the bifurcation 
parameter A is lower (left) or higher (right) than the critical value A c . This bifurcation parameter can be 
either the density or the inactivation variable of a cooperative channel. Other adaptation variables are set to 
constant values chosen in physiological ranges (see text). For A lower (resp. higher) than the critical value A c , 
all models exhibit electrophysiological signatures of competitive (resp. cooperative) excitability. Numerical 
values of the parameter A in the different plots are as follows. Thalamic relay cell: left hc a ,T = 0, right 
h>Ca,T = 0.2. Thalamic reticular cell: left hc a ,T = 0, right hc a ,T = 0.4. Aplysia R15 neuron: left 
§Ca.L = 10~ 6 mS/cm 2 , right gca,L = 0.00AmS/cm 2 . GC neuron: left hc a ,N = 0, right hc a ,N — 0.3. 
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Figure 10. Thalamocortical reticular and relay cell signaling relies on dynamical switches 
between competitive and cooperative excitability, through the regulation of T-type calcium 
channel inactivation. A. Electrophysiological responses of a full TC reticular cell model [5] to a step input 
of excitatory current starting from a hyperpolarized (left) or depolarized (right) steady-state. An initial 
hyperpolarization of the cell releases T-type calcium channel inactivation, inducing a switch from competitive 
to cooperative excitability responsible for the two distinctly different electrophysiological responses. Feedback 
interconnections illustrate the switch between competitive and cooperative excitability induced by the 
regulation of T-type calcium channel inactivation (the darker the arrows, the stronger the feedback loop). B. 
Electrophysiological responses of the same TC reticular cell model to a step input of excitatory current of 
higher amplitude, starting from a hyperpolarized resting state. The figure shows variations of membrane 
potential (top) and T-type calcium channel inactivation (bottom) over time. The neuron dynamically switches 
from cooperative to competitive excitability during the excitation, due to T-type calcium channel inactivation. 
This results in the generation of a transient burst followed by low-frequency single-spike firing. 
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Figure 11. The heterogeneity of dopaminergic (DA) neuron L-type calcium channel density is 
responsible for an heterogeneity in their firing pattern. A. Extracellular recordings of substantia 
nigra of two DA neurons from rat brain slices in control condition (top), during SK channel blockade (middle), 
and during concomitant blockade of SK channel and L-type calcium channel (bottom). B. Coefficient of 
variations of neuron firing patterns in each condition. Black point represents values for each neuron, bars 
represents the mean values and the standard deviation. C. Feedback interconnections illustrating the 
heterogeneity in DA neuron firing pattern after SK channel blockade and in the effect of the L-type channel 
blocker nifedipine (the darker the arrow, the stronger the feedback). Although DA neurons exhibit similar 
firing patterns in control conditions, SK channel blockade uncovers an heterogeneity, some DA neuron 
exhibiting electrophysiological signature of competitive excitability, namely single-spike firing, others of 
cooperative excitability, namely burst firing. A L-type calcium channel blockade solely affect the firing pattern 
of bursting neurons, inducing a switch from cooperative to competitive excitability, as predicted above. 
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Supplementary material 



Construction of the transcritical bifurcation in generic conductance- 
based models 

Suppose that, given a bifurcation parameter fi e {cji, Ei}i e x, there exists a a solution (V, /i) = (Vtc, Mtc) 
satisfying the two degeneracy conditions (7) and (8). Then we claim that, posing 

Itc = Y,^ m 1UVTc)h^(VTc)(V T c - E t ), 

i 

which solves the fixed point equation, the model is at a transcritical bifurcation. The existence of the couple 
(Vtc, fJ-Tc) in the physiological range is shown for specific examples in the main document. 

The center space E c associated to the degeneracy conditions (7) and (8) is spanned by the vector 



(l,k), k 



dx Q 



dV 



v=v T 



where x runs all fast and recover gating variables in the model, that is k is the vector whose elements are the 
slopes of the (in) activation functions of all the fast and recovery gating variables calculated at V = Vtc- The 
associated center manifold M c is exponentially attractive. Indeed, it can easily been shown that, when (7) and 
(8) are satisfied, the remaining nonzero eigenvalues of the Jacobian of (10) are all negative. To the first order, 
the dynamics on M. c is given by 

C m V c = -J29i(m i ,o {V T c) + k mi (V c - V TC )) ai (h i>00 {V TC ) + k hi {V c - V TC )) b \V c - E t ) + I app 
iei 

Consider the affine reparametrization 



M Mj I app — I app 



dV 



(P - fJ-Tc)- 



(S16) 



M=Mtc 
V=V TC 



It is easy to show that, in the affine reparametrization (£16 1, the center manifold dynamics satisfy 



dVc 
dV 



V=Vtc 



dVc 
dp, 



fJ-=fJ-TC 

V=-Vtc 



corresponding to the defining condition of a transcritical bifurcation (see [U Page 367]). 
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